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BACKGROUND 

During  the  Gemini,  Apollo,  and  Sky lab  programs,  cryogenic  systems  were 
used  as  the  source  for  fuel  cell  energy  production.     Throughout  these 
programs  the  stored  cryogens  became  one  of  the  main  limiting  consumables 
involving  vehicle  lifetime.     The  critical  nature  of  the  cryogenic  system 
required  very  accurate  premission  analyses.     Not  only  were  quantity 
calculations  required  but  pressure  range  evaluation  was  also  made  to  ensure 
that  mission  operations  could  be  safely  performed. 

A  potential  major  improvement  in  Space  Shuttle  operations  involves  the 
use  of  slush  cryogens.     Slush  is  the  term  used  to  describe  a  mixture  of  solid 
and  liquid  phases  of  a  cryogen.     The  advantage  of  slush  over  liquid  cryogens 
is  a  higher  stored  mass  for  the  same  occupied  volume  storage  tank.  Also, 
the  latent  heat  of  the  solid  in  the  slush  increases  storage  time  before 
venting  occurs.     This  improved  capacity  and  storage  time  make  possible 
increased  space  flight  time  and  distance. 

The  Cryogenics  Division  of  the  National  Bureau  of  Standards  provided 
the  ?VT  data  and  the  equations  of  state  which  were  used  during  the  Apollo  and 
Sky lab  programs  under  the  sponsorship  of  NASA.     The  basic  data  and  computer 
programs  are  extended  in  the  present  report  to  provide  the  support  necessary 
to  an  increased  effort  in  slush  research  and  mission  analysis.     This  software 
technology  is  now  available  for  all  future  space  programs  and  is  expected  to 
have  industrial  application.  ... 

Paralleling  the  basic  data  effort  NASA  has  been  a  prime  driver  in 
funding  the  development  of  slush  cryogen  technology.     A  slush  hydrogen 
program  conducted  from  1965  to  1970  at  NBS  investigated  several  production 
methods,  measured  the  pressure  drop  in  flowing  slush  hydrogen,  determined 
maximum  solid  fraction  and  investigated  several  instruments  to  measure  solid 
fraction  .     An  experimental  program  currently  underway  at  NBS  is  expected  to 
demonstrate  the  safe  and  reliable  production,   storage  and  handling  of  slush 
hydrogen  and  oxygen  in  ground  support  equipment  at  NASA  launch  facilities. 


* 

See  page  21,   22  of  Nat.   Bur.   Stand.   Technical  Note  664  for  references. 
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SYMBOLS 


P 

pressure 

dP 
dT 

derivative  of  vapor  pressure  or 

melting  pressure 

T 

temperature 

T' 

temperature  shifted  to  a  new  scale 

P 

density 

s 

specific  heat  at  constant  pressure 

V 

volume  - 

1/p  •■ 

specific  heat  at  constant  volume 

H 

enthalpy 

AH 

enthalpy  difference 

S 

entropy 

AS 

entropy  difference 

U 

internal 

energy 

c 

SAT 

specific  heat  of  saturated  solid  or 

saturated  liquid 

Q 

quality 

A, 

B, 

C ,   a ,   c  , 

P 

o 

coefficients 
equations 

in 

vapor  pressure  or  melting  pressure 

SUBSCRIPTS 


liquid 
vapor 
solid 
t 

melt 


liquid  phase 
vapor  phase 
solid  phase 
triple  point 
melting  line 


fusion 

sublimation 

vaporization 


UNITS 

The  primary  variables  in  the  computer  programs  are 

Pressure  in  atmospheres 

Density  in  moles/liter,  and 

Temperature  in  Kelvin. 
Conversions  to  other  SI  units  and  units  normally  used  in  applied  problems 
are  given  in  Appendices  E  and  F. 
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THE  THERMODYNAMIC  PROPERTIES  OF  SLUSH  HYDROGEN  AND  OXYGEN 


Hans  M.  Roder 

Cryogenics  Division 
Institute  for  Basic  Standards 
National  Bureau  of  Standards 
Boulder,   Colorado  80302 


The  thermodynamic  properties  of  hydrogen  and  oxygen  have  been  calculated 
for  temperatures  both  below  and  above  the  triple  point.     Values  of  the  various 
physical  properties  required  for  the  computations  are  either  taken  from  the 
literature,  or  are  extrapolated.     If  extrapolated,  the  extrapolations  are 
based  on  the  known  behavior  of  other  simple  fluids.     The  results  are 
presented  in  the  form  of  computer  programs  which  cover  two  regions  of  the 
phase  diagram,   solid-vapor,   and  solid-liquid.     Input  to  the  programs  is 
temperature  and  quality.     The  properties  returned  include  pressure,  density, 
enthalpy,  entropy  and  internal  energy.     The  present  programs  for  slush 
hydrogen  and  slush  oxygen  are  an  extension  of  and  depend  upon  the  equations 
of  state  developed  previously  for  these  gases. 

Key  words:     Computer  programs;   enthalpy;   entropy;   hydrogen;   internal  energy; 
liquid;   oxygen;   PVT;   quality;    slush;    solid;  vapor. 

1 .  INTRODUCTION 

Slush  is  the  term  used  to  designate  a  mixture  of  solid  and  liquid 
phases  of  the  gas  under  consideration.     There  are  two  major  reasons  for  the 
current  interest  in  slush.     First,   since  solid  is  more  dense  than  liquid, 
slush  offers  a  way  to  carry  a  larger  quantity  of  the  cryogen  in  a  container 
of  fixed  volume,   as  compared  to  the  liquid  alone.     Second,   since  the  solid 
melts,  absorbing  heat,   slush  offers  a  way  to  lengthen  the  time  that  the 
cryogen  can  be  stored  in  a  fixed  volume  container  before  venting  has  to 
occur. 

We  have  included  the  solid-vapor  region  in  this  study  because  slush  is 
currently  produced  by  a  freeze-thaw  process.     This  means  the  properties  of  the 
vapor  at  temperatures  just  below  the  triple  point  are  also  of  concern.  Further, 
the  properties  of  this  two  phase  region  are  of  interest  in  another  application, 
subliming  refrigerators,  which  are  used  to  cool  infrared  detectors. 

The  properties  of  concern  in  this  study  are  temperature,  pressure, 
density,   enthalpy,   entropy  and  internal  energy.     Internal  energy  was  added 
to  the  set  because  the  principal  applications  at  present  involve  containers 
of  fixed  volume.     The  regions  of  the  phase  diagram  covered  in  this  report, 
solid-vapor  and  solid-liquid,   are  both  two  phase  regions.      In  analogy  -o  the 
familiar  quality  used  in  liquid-vapor  problems,   we  have  included  qualicy  as 
a  parameter  in  this  study. 


The  thermodynamic  properties  described  in  this  report  are  consistent 
with  and  should  be  viewed  as  an  extension  of  the  two  equations  of  state 
developed  previously,  the  equation  of  state  for  hydrogen  by  Roder  and 
McCarty   (1975) ,   and  the  equation  of  state  for  oxygen  by  Weber   (1977a,   1977b) . 
The  single  phase  regions  of  each  gas  for  liquid,  and  fluid  gas  or  vapor  as 
well  as  the  liquid-vapor  two  phase  region  are  defined  by  these  equations  of 
state.     The  present  slush  routines  call  on  the  equations  of  state  to  define 
liquid  properties  along  the  melting  line  and  vapor  properties  on  the  vapor 
boundary  of  the  solid-vapor  region. 

2.     BASIC  THERMODYNAMIC  CONSIDERATIONS 

2.1  General 

Consider  a  pressure  temperature  phase  diagram  shown  schematically  in 
figure  1.     We  recognize  the  melting  line,   the  vapor  pressure  curve  of  the 
liquid  and  the  vapor  pressure  of  the  solid.     These  curves  define  the 
boundaries  of  the  various  phases,   solid,   liquid,   and  vapor.     They  intersect 
at  the  triple  point.     It  is  from  a  plot  such  as  this  that  the  triple  point 
takes  its  name.     What  this  phase  diagram  does  not  show  is  the  variation  of 
density.     One  unique  aspect  of  any  two  phase  regime  liquid-vapor  solid-vapor 
or  solid-liquid  is  that  the  pair     P,  T    on  a  phase  boundary  encompasses  a 
large  range  of  densities.     We  could  consider  a  plot  in  the  three  variables 
pressure,  density  and  temperature,  however,  consider  instead  the  more 
familiar  T-S  diagram  shown  schematically  in  figure  2.     This  diagram  is  also 
a  phase  diagram.     We  recognize  the  single  phase  regions  solid,  liquid  and 
vapor.     We  can  also  see  the  various  two  phase  regions  where  two  phases  are 
in  equilibrium.     In  a  diagram  of  this  kind  the  triple  point  gets  stretched 
into  a  line.     The  density  variation  can  also  be  included.     For  example,  the 
density  at  point  1  in  figure  2  corresponds  to  the  density  of  the  solid  at  the 
triple  point,  point  2  to  the  liquid  triple  point  density,  and  point  3  to  the 
vapor  density  at  the  triple  point. 

The  regions  already  defined  by  the  major  equations  of  state  are  shown 
shaded-     We  will,   in  this  report,  cover  the  two  phase  regimes  solid-liquid 
and  solid-vapor.     If  we  assume  that  the  P-T  relation  of  the  boundary  is 
available,   then  to  define  a  region  we  need  only  define  the  properties  at  the 
boundaries.     For  the  true  slush  region,   solid-liquid,  we  would  have  to  know 
density,  enthalpy,   entropy  and  internal  energy  at  both  the  liquid  and  solid 
boundaries.     The  liquid  side  is  already  defined  for  us  by  the  intersection 
of  the  equation  of  state  and  the  melting  pressure  curve.     It  remains  to 
define  the  other  boundary  in  this  case  the  solid  side.     The  actual  approach 
depends  on  what  pieces  of  experimental  data  are  available  to  guide  us,  and 
these  are  discussed  in  detail  in  sections  that  follow. 
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Figure  1.     Pressure  vs.  Temperature,  schematic. 
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A  fundamental  relation  governs  the  various  two  phase  regimes.  This 
equation  is  called  the  Clausius-Clapeyron  equation,   and  we  will  use  it 
extensively.     Written  out  for  each  of  the  three  two  phase  regimes  the 
equation  reads  as  follows: 

Solid-Liquid  Solid-Vapor  Liquid-Vapor 

dP  _  AH   (Fusion)   dP  _  AH   (Sublimation)       dP  _  AH  (Vaporization) 

dT       T(V,  .      .-,-V     ,  .  ,)      dT       T(V  -V     ,.,)        dT       T  (V  -V,.  .,) 

^  liquid     solid  vapor     solid  vapor  liquid 

The  equation  relates  a  latent  heat,  the  pressure-temperature  derivative 
at  the  appropriate  phase  boundary  and  the  volume  change  corresponding  to  the 
phase  change  in  question.     If  the  available  data  are  redundant,   then  we  are 
in  a  position  to  make  thennodynamic  consistency  checks.     Quite  often  the  list 
of  available  data  is  minimal  and  the  opportunity  for  consistency  checks  is 
lost.     When  elements  of  data  are  missing  we  are  forced  to  calculate  the 
unknown  pieces  of  information  from  the  Clausius-Clapeyron  equation.  For 
example,  the  volume  change  on  fusion  is  usually  not  measured  experimentally. 
The  heat  of  fusion  is  normally  measured  only  at  the  triple  point,  not  at 
elevated  temperatures  or  pressures.     If  we  wish  to  know  the  volume  change  on 
fusion,  say  at  the  triple  point,   then  we  must  calculate  it  from  the 
experimentally  measured  heat  of  fusion  and  the  derivative  of  the  melting 
pressure  equation. 

2.2     The  Quality,  Q 

The  concept  of  quality  is  familiar  to  most  engineers  who  have  solved 
liquid-vapor  problems.     The  conventional  quality  is  defined  as  the  percent 
(by  mass)   of  the  vapor  phase  present,  or: 

.,)   for  volume,  and  similarly 
liquid  ^ 

liquid^ 

,  .      .  , )  and 
liquid 

liquid 

Since  it  is  easy  to  make  an  error  we  caution  that  the  expression  in  teims  of 
density  is : 


Q  = 

(V 

-V,  .      •  .)/(V 

liquid  vapor 

-V 

Q  = 

(H 

-H,  .      •  -,)/  (H 

liquid  vapor 

-H 

Q  = 

(S 

-S,  .      •  -,)/(S 

liquid  vapor 

-S 

Q  = 

(U 

-u,  .     •  -,)/ (U 

liquid  '  vapor 

-u 

Q  = 


Pi  •      ■ J  •  P 

liquid  vapor 

p  ' vapor 


/  (  P 1  •      ■   -  P  ) 
liquid  vapor 


From  these  equations  it  is  apparent  that  Q  =  0.0  on  the  liquid  side,  i.e., 
the  denser  side,   and  that  Q  =  1.0  on  the  vapor  side,   i.e.,   the  less  dense 
side. 
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For  the  other  two  phase  regimes  solid-liquid  and  solid-vapor  we  will 
use  a  defined  quality     Q     defined  in  a  manner  quite  analogous  to  that  for 
liquid-vapor  above-     Q    will  be  0.0  on  the  dense  side,   i.e.,   the  solid,  and 
Q    will  be  1.0  on  the  less  dense  side.     The  equations  read  as  follows: 

Solid-Liquid  Solid-Vapor 

Q    =     (V-^'solid^/^^liquid-^solid^       Q    =  (^-^solid^/^^vapor-^solid^ 

Q     =      («-«S0lid)/(«liquid-«S0lid)        Q     =  (H-H3^,.^)/(H^^p^^-H3^^.^) 

Q     =      (S-S     T.,)/(S,  .      .,-S     ,.,)        Q     =      (S-S     T.,)/(S  -S     ,  .  ,) 

^  solid  '      liquid     solid  solid  '     vapor  solid 

Q     =      (U-U     ..,)/(U,  .      .-,-U  Q     =      (U-U     ..,)/(U  -U     ,  .  ,) 

solid^      liquid     solid  solid        vapor  solid 


3.     THE  SOLID-LIQUID  REGION  FOR  HYDROGEN 

Prior  work  on  the  thermodynamic  properties  of  this  region  include  the 
experimental  work  and  the  calculations  by  Dwyer,  et  al.    ( 1965a, b)   and  by 
Cook,  et  al.    (1965) .     In  the  present  report  the  thermodynamic  properties 
along  the  liquid  boundary  of  this  region  are  defined  from  the  intersection 
of  the  32-term  modified  Benedict-Webb-Rubin   (MBWR)   equation  of  state  as 
presented  by  Roder  and  McCarty   (1975)   with  the  melting  pressure  equation 
(Goodwin,   1962) .     The  elements  of  data  available  to  describe  the  region  are 
limited.     In  the  sections  that  follow  we  describe  the  reference  values  chosen 
for  the  property  calculations,  the  melting  pressure  equation  of  the  solid, 
the  heats  of  fusion  at  elevated  pressures,  the  property  calculations,  and 
thermodynamic  consistency  checks. 

3.1    The  Reference  Values 

Since  the  available  set  of  data  is  limited  we  are  forced  to  calculate 
changes  in  properties  for  this  region  from  the  Clausius-Clapeyron  equation. 
To  obtain  values  for  the  solid  boundaries  these  changes  are  applied  from 
the  liquid  boundary.     In  effect,   for  a  given  temperature  the  values  obtained 
from  the  32-term  equation  of  state  and  the  melting  line  are  used  as  starting 
point,  or  point  of  reference.     At  the  triple  point  temperature  the  procedure 
leads  to  values  for  the  solid  that  are  used  as  reference  point  for  the 
solid-vapor  region   (see  section  4.1). 

3.2     The  iMelting  Pressure  Equation 

The  equation  used  is  the  one  developed  by  Goodwin   (1962),   except  that 
we  shift  the  temperature  scale  by  0.003  K  to  match  the  triple  point 
temperature  and  pressure  of  the  liquid  vapor  pressure  curve  used  with  the 
32-term  MBWR.     The  shift  insures  that  the  melting  pressure  equation  and  the 
liquid  vapor  pressure  equation  intersect  at  the  triple  point,   13.800  K  and 
0.0695  atm,   which  are  the  values  chosen  for  the  basic  32-term  equation  of 
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state.     With  the  temperature  shift  Goodwin's  equation  becomes: 
P  -  P 

T^-^    =     A  e  -  ^  ^  BT' 

where  t'   =  T  +  0.003  in  kelvin,  A  =  30.3312,   a  =  5.693,  B  =  2/3,  P^  =  0.0695 
atm  and        =  13.803. 

3.3     The  Heats  of  Fusion  at  Elevated  Pressures 

These  values  were  determined  experimentally  by  Dwyer,  et  al.  (1965b). 
The  authors  represent  the  experimental  results  by  a  straight  line  in  pressure. 
We  select  the  value  of  117.277  J/mol    (28.03  cal/mol,  Woolley,  et  al.,  1948) 
for  the  heat  of  fusion  at  the  triple  point,  partially  to  remain  consistent 
with  the  earlier  thermodynamic  calculations,  and  partially  because  this  value 
is  quite  close  to  the  average  of  all  experimental  determinations   (see  table 
9-3,   Roder,   et  al.,   1973).     Shifting  the  Dwyer,   et  al.   equation  to  pass 
through  the  selected  value  the  equation  used  for  the  heats  of  fusion 
becomes : 

AHp     .         =      (28.02693164  +  0.044149P)    •  4.184 
fusxon 

where  AH^     •       is  in  J/mol  and     P     in  atm.     The  uncertainty  in  AH^; 

fusion  '  ^  fusion 

remains  essentially  as  estimated  by  Dwyer,   et  al.    (1965b),   around  1%. 

3.4     Property  Calculations 

To  find  the  values  of  the  various  properties  on  the  solid  boundary  we 
calculate  the  volume  and  entropy  changes  from  the  Clausius-Clapeyron  equation. 
The  melting  pressure  equation  and  the  heats  of  fusion  were  given  in  the 
previous  sections.     The  relations  are: 

P     =     P,,elt(T^  . 
dF  ^^^melt^ 


dT  dT 


AH  (fusion) 


V  =     V,.      --.-AV^      .       with  AV^      .  -   

solid  liquid       fusion  fusion  dP 

dT 


H     T  ■  J     ~     Ht  .      •      —  AHj- 
solid  liquid  fusion 


,  .  -,     =     S,  .      .  -,   -  AS^      .      ,   with  AS^      .  =     AH^      .  /T. 

solid  liquid  fusion  fusion  fusion 
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The  internal  energy  for  the  solid  follows  from  the  thermodynamic  relation 
U  =  K  -  PV. 

3.5     Thermodynamic  Consistency  Checks 

Since  the  available  set  of  data  is  limited,   the  checks  that  can  be  made 
are  limited  also.     We  can  examine  the  change  of  volume  on  fusion  as  a 
function  of  temperature.     We  expect  the  volume  change  to  decrease  as  the 
temperature  increases,   it  does.     The  values  obtained  from  the  heats  of  fusion 
through  the  Clausius-Clapeyron  equation  are  in  excellent  agreement  with  the 
experimental  volume  changes  published  by  Bartholome   (1936)   for  normal 
hydrogen  as  shown  in  table  1. 


Table  1.     Volume  Changes  on  Fusion  at  Elevated  Temperatures 

T  Bartholome  This  Report  % 

13.96  2.83  2.81  .58 

16.43  2.31  2.26  2.15 

18.24  2.03  2.02  .61 


Considering  nitrogen  and  argon  we  would  expect  that  the  entropy  change 
on  fusion  for  hydrogen  should  decrease  as  the  temperature  increases. 
Quantitative  comparisons  cannot  be  made;  however,  this  expectation  is  at 
least  qualitatively  fulfilled. 

4.      THE  SOLID-VAPOR  REGION  FOR  HYDROGEN 

Prior  work  on  the  thermodynamic  properties  of  this  region  include 
calculations  by  Mullins,   et  al.    (1961),   the  T-S  chart  by  Sindt  and  Mann  (1966) 
and  the  survey  by  Roder,  et  al.    (1973) .     The  thermodynamic  properties  along 
the  vapor  boundary  of  this  region  are  available  from  the  intersection  of  the 
32-term  MBWR   (Roder  and  McCarty,   1975)   and  the  solid  vapor  pressure  curve 
defining  the  boundary  of  the  region.     The  elements  required  to  complete  the 
description  of  the  region  are  the  reference  point  for  all  properties,  the 
vapor  pressure  equation  of  the  solid,   the  density  of  the  solid,  and  the  heat 
capacity  of  the  solid. 

4.1     The  Reference  Point 

we  wish  to  integrate  this  region  in  a  thermodynamically  consistent  way 
with  the  basic  32-term  equation  of  state.     The  simplest  thing  we  can  do  is  to 
refer  the  properties  of  the  solid  at  the  triple  point  to  those  obtained  from 
the  32-term  MBWR  for  the  liquid  at  the  triple  point.     Specifically,  we 
require  that  at  the  triple  point  temperature: 
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vapor  pressure  of  the  solid    =    vapor  pressure  of  the  liquid 
solid  Ixquid  fusion 

solid  liquia  fusion 

S  =     S  —  AS 

solid  liquid  fusion 

Internal  energy  of  the  solid  follows  through  the  thermodynamic  relation 

U  =  H  -  PV.     AHc     .      ,   the  heat  of  fusion  selected  for  the  triple  point  is 
fusion  ^ 

117.277  J/mol   (see  section  3.3).     AS^  ,  the  entropy  change  on  fusion  is 

fusion  '^-^  ^ 

defined  thermodynamically  as  AH^  /T.     AV_     .      ,   the  volume  change  on 

^  fusion'  fusion  ^ 

fusion  is  calculated  from  the  Clausius-Clapeyron  equation  and  the  melting 

pressure  derivative   (see  section  3.4).     The  result  for  the  triple  point  with 

=  13.  800  K  is  AV^  =2.86.  cm'^/mol. 

t  fusion  4 

Having  established  a  reference  point  or  benchmark  we  can  look  next  at 

the  description  of  the  remainder  of  the  region. 

4.2     Vapor  Pressure  of  the  Solid 

One  of  the  drawbacks  of  the  earlier  calculation  by  Mullins,  et  al.  (1961) 
is  that  the  authors  presented  no  explicit  equation  for  the  vapor  pressure  of 
the  solid,   indeed  they  could  not  because  one  can  either  force  exact 
thermodynamic  consistency  or  use  an  explicit  vapor  pressure  equation,  but  not 
both.     Experimental  measurements  of  the  vapor  pressure  of  solid  hydrogen 
exist  to  quite  low  temperatures   (Woolley,   et  al.,   1948;   Borovik,   et  al . ,  1960; 
and  Harrison,   et  al.,   1962).     However,  we  will  use  these  values  only  as  a 
check  on  the  extrapolation  of  the  vapor  pressure  curve  to  low  temperatures. 

The  form  of  the  vapor  pressure  equation  is  dictated  by  thermodynamics. 
Integration  of  the  Clausius-Clapeyron  equation  assuming  the  vapor  phase  to 
be  an  ideal  gas  leads  to  the  Kirchhof f -Rankine  equation: 

log  P     =     A/T  +  B  log  T  +  C 

We  use  this  equation  but  subject  it  to  the  following  constraints: 

1.  P^  =  0.0695  atm,   T^  -  13.800  K.     To  insure  that  this  equation 
intersects  the  melting  equation  and  the  liquid  vapor  pressure  equation 
at  the  triple  point  which  is  defined  by  the  liquid  vapor  pressure 
equation  used  with  the  32-term  MBWR. 

dP 

2.  ^    at  the  triple  point  temperature  =  0.04566744  atm/K.     This  constraint 
defines  the  derivative  dP/dT  of  the  solid  vapor  pressure  curve  at  the 
triple  point  temperature,   and  through  the  Clausius-Clapeyron  equation 
the  heat  of  sublimation.     The  value  of  the  derivative  is  chosen  so  that: 

AH^  ,  ,  .      ^.  =     AH^  +     AH  .  ^. 

.       .  sublimation  fusion  vaporization 


where  IH-  is  the  value  chosen  in  sections  3.3  and  4,1,   and  the 

rusion 

IH  is  determined  by  the  32-term  MBWR  and  its  companion 

vaporizarion  ^ 

liquid  vapor  pressure  equation. 

3.       That  the  vapor  pressure  at  10  K  be  0.00255940  atm.     The  last  requirement 
is  the  best  compromise  to  insure  thermodynamic  consistency  of  the  heat 
of  sublimation  down  to  10  K.     The  point  is  discussed  in  more  detail  in 
the  section  on  thermodynamic  consistency  below. 

With  3  degrees  of  freedom  and  three  constraints  the  coefficients  of  the 
solid  vapor  pressure  equation  are  uniquely  determined.     They  are: 

A    =     -90.77568949         B     =     2.489830940         and     C     =  4.009857354 

where    P     is  in  mm  Hg  and    T    is  in  K. 

Intercomparison  to  the  available  experimental  vapor  pressures  shows  the 
maximum  deviation  in  pressure  down  to  10  K  is  less  than  1%.     The  simple 
equation  cannot,  however,  represent  the  data  down  to  low  temperatures,  say 
3  or  4  K.     In  this  temperature  range  the  departure  between  predicted  and 
experimental  vapor  pressures  increases  to  nearly  200%.     Additional  terms  in 
the  vapor  pressure  equation  improve  the  representation  of  the  vapor  pressures 
but  not  the  heats  of  sublimation  calculated  from  the  vapor  pressure  derivative. 

4.3     Density  of  the  Solid 

From  the  survey  by  Roder ,   et  al.    (197  3)   we  select  the  solid  volumes 
published  by  Ahlers    (1963)   primarily  because  we  expect  these  values  to  be  in 
accord  with  the  heat  capacity  measurements  selected  in  the  next  section.  More 
recent  publications    (Anderson  and  Swenson,   1974)   prefer  values  which  differ 
by  as  much  as  1.5%  from  the  ones  selected  especially  at  the  lower  tempera- 
tures.    Ahlers'   results  can  be  represented  by  a  parabolic  function  from  6  to 
13.3  K.     We  shift  this  function  to  pass  through  the  designated  value  of  the 
solid  density  at  the  triple  point.     The  resulting  relation  is: 

rr^   =     V  =     22.  77963478   -  0.  031T     +  0.005T^ 

^solid  ^""^'-^ 

with  V     , . ,  in  cm"/mol  and     T     in  K.     We  estimate  that  the  values  of  the 
sclxa 

solid  should  be  in  error  by  no  more  than  0.4%  at  the  triple  point,  and  that 
the  error  could  increase  -o  as  ~uch  as  1.5%  at  4  K  or  lower. 

4,4     Heat  Capac:ity  of  -che  Solid 

We  require  heat  capacity  values  for  the  solid  at  saturarion .     The  values 
selected  were  derived  by  Ahlers    (1963)    from  his  measurements  of  C^.  in  -he 
3::.ngle  phase  solid.     We  orefsr  these  values  to  earlier  measuremen-cs  because 
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a  closed  loop  third  law  check  of  entropy  using  these  values  is  much  more 
consistent  than  checks  using  the  earlier  measurements   (See  Roder,  et  al., 
1965) .     Since  several  conventional  representations  do  not  seem  to  apply  to 
the  values  for  hydrogen,  we  settled  on  a  simple  interpolating  polynomial 
which  returns  Ahlers'  values  to  within  0.3%  from  14  K  down  to  5  K»  Maximum 
deviations  down  to  3  K  are  1.65%,  but  the  polynomial  should  not  be  used  below 
3  K.     The  power  series  used  to  represent  C^^^  is: 

C„-^     =     -101.     +     87.75  T  -  26.  95833333         +  4.25  T"^  -  0.041666667 

with  Cg^^  in  milli  Joules/mole  and     T     in  K. 

We  use  the  heat  capacity  in  the  thermodynamic  relations: 

^«    =    ^l""  Sat  =    /^^  dT 

T 

to  calculate  the  changes  in  enthalpy  and  entropy  from  the  reference  values 
for  the  solid  at  the  triple  point. 

4.5    Thermodynamic  Consistency  Checks 

Once  the  vapor  pressure  of  the  solid  is  defined,  we  can  perform  a 
consistency  check  in  the  following  way.     We  calculate  the  heat  of  sublimation 
first  by  taking  the  difference  between  vapor  and  solid  enthalpies  at  the 
boundary,   and  second  by  calculating  it  from  the  Clausius-Clapeyron  equation 
using  solid  and  vapor  volumes  and  the  derivative  of  the  solid  vapor  pressure 
curve.     The  enthalpy  of  the  vapor  is  obtained  from  the  32-term  equation  of 
state,  we  note  that  these  enthalpies  depend  on  the  specific  heat  of  the  ideal 
gas  which  is  not  defined  below  10  K  in  the  equation  of  state  programs.  The 
enthalpy  of  the  solid  is  obtained  from  the  integration  of  the  heat  capacity 
of  the  solid.     The  difference  of  these  two  values,   the  heat  of  sublimation, 
calculated  in  this  way  is  not  particularly  sensitive  to  the  value  of  pressure 
chosen.     In  the  other  calculation,   however,   the  heat  of  sublimation  depends 
on  the  vapor  volume  obtained  from  the  equation  of  state  and  on  the  derivative 
of  the  solid  vapor  pressure  curve  and  is  quite  sensitive  to  changes  in 
pressure.     We  select  the  value  of  the  pressure  at  10  K,   the  third  constraint, 
to  minimize  the  departures  between  the  heats  of  sublimation  calculated  by  the 
two  different  methods  between  13.8  and  10  K.     The  intercomparison  of  the 
heats  of  sublimation  is  shown  in  table  2;   the  departures  do  not  exceed  0.1%. 
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Table  2.     Intercomparison  of  the  Heats  of  Sublimation  of  Hydrogen 


T  AH   (sublimation)    from        AH   (sublimation)  from  the  Difference 

K  the  boundaries  J/mol          Clausius-Clapeyron  Eq.  in  percent 

J/mol  J/mol 

10  961.234  960.300  0.10 

11  979.262  979.169  0.01 

12  996.088  996.584  -0.04 

13  1011.351  1011.900  -0.05 

13.8                      1022.137  1022.137  defined  to  be  the  same 


5.      THE  SOLID-LIQUID  REGION  FOR  OXYGEN 

With  the  exception  of  the  melting  curve  and  the  heat  of  fusion  at  the 
triple  point  prior  work  on  the  properties  of  this  region  is  non-existent. 
We  are  therefore  forced  to  estimate  all  of  the  thermodynamic  properties  from 
the  Clausius-Clapeyron  equation  using  other  simple  gases  as  model  fluids. 
The  thermodynamic  properties  along  the  liquid  boundary  of  this  region  are 
defined  from  the  intersection  of  the  32-term  MBWR  equation  of  state  and  the 
melting  pressure  equation  presented  by  Weber   (1977b).     Weber's  version  of 
the  32-term  MBWR  rather  than  an  earlier  version   (Stewart,  et  al.,   1972)  is 
used  because  it  includes  new  values  of  PVT   (Weber,   1977a)   which  go  to  high 
pressures   {"^  800  bar).     An  additional  advantage  of  Weber's  formulation  is 
that  the  vapor  pressure  curve  and  the  melting  line  intersect  at  the  defined 
triple  point.     In  the  sections  that  follow  we  describe  the  reference  values 
chosen  for  the  property  calculations,   the  melting  pressure  equation  of  the 
solid,   and  the  property  estimation. 

5.1     The  Reference  Values 

Since  the  available  set  of  data  is  limited,  we  are  forced  to  calculate 
changes  in  properties  for  this  region  from  the  Clausius-Clapeyron  equation. 
To  obtain  values  for  the  solid  boundaries  these  changes  are  applied  from  the 
liquid  boundary.     In  effect,   for  a  given  temperature  the  values  obtained 
from  the  32-term  equation  of  state  and  the  melting  line  are  used  as  a  starting 
point,   or  point  of  reference.     At  the  triple  point  temperature  the  procedure 
leads  to  values  for  the  solid  that  are  used  as  a  reference  point  for  the 
solid-vapor  region   (see  section  6.1).     The  enthalpy  departure  at  the  triple 
point  has  been  estimated  by  Weber   (1977b)    to  be  about  3  J/mcl. 

5.2     The  iMelting  Pressure  Equation 

The  equation  used  is  the  one  developed  by  Weber   (1970)   as  modified  by 
him   (Weber,   1977a, b)    to  pass  through  the  defined  triple  point.     The  equation 
reads 
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with        =  54.359  K,         =  0.0014451  atm,   P^  =  2637.33  atm,  and  c  =  1.769. 

5.3     Property  Estimation 

To  estimate  values  of  the  various  properties  on  the  solid  boundary  of 
this  region  we  will  use  the  Clausius-Clapeyron  equation.     To  do  this  we  must 
make  an  educated  guess  as  to  how  the  heat  of  fusion,  or  the  volume  change  on 
fusion,  or  the  entropy  change  on  fusion  depends  on  pressure   (or  temperature) . 
Considering  the  known  behavior  of  hydrogen,   argon,   and  nitrogen  we  postulate 
that  oxygen  should  obey  the  following  qualitative  criteria: 

1.  The  heat  of  fusion  should  increase  as  either  pressure  or  temperature  is 
increased, 

2.  The  entropy  of  fusion  should  decrease  as  either  pressure  or  temperature 
is  increased,  and 

3.  The  volume  change  on  fusion  should  decrease  as  either  pressure  or 
temperature  is  increased. 

Direct  scaling  of  the  volume  change  for  oxygen  from  hydrogen,  argon,  or 
nitrogen   (see  page  53,   Sindt,   Ludtke  and  Roder ,   1970)    is  not  successful;  one 
of  the  criteria  remains  unfulfilled.     Adopting  a  sort  of  an  average  between 
the  two  extreme  assumptions 

1.       AH^     .        (T  or  P)   =  constant,  and 
fusion   ^  '  ' 

^*       '^^fusion   ^"^  °^        ~  constant 

allows  us  to  fulfill  all  three  criteria  over  a  range  of  experimentally 

plausible  pressures,   say  up  to  5000  psia.     As  for  hydrogen  we  express  the 

heat  of  fusion  of  oxygen  as  a  linear  function  of  pressure 

^^fusion     "      (106.3  +  0.01126  P)    •  4.184 

where  ^^lf^sj_Qj^        iri  J/mol  and  P  is  in  atm.     The  value  106.3  cal/mol  is 
taken  at  the  triple  point  from  the  paper  by  Giauque  and  Johnston   (1929) . 
The  uncertainty  in  the  estimated  heat  of  fusion  should  be  no  larger  than  the 
variation  estimated  for  it,   i.e.,   a  maximum  of  6%  at  the  highest  pressures. 

With  an  estimate  for  the  heat  of  fusion  in  hand  the  property  calculations 
follow  the  relations  given  in  section  3.4.     Obviously,   consistency  checks 
cannot  be  made  for  this  region. 

■  ■  '    6.      THE  SOLID-VAPOR  REGION  FOR  OXYGEN 

Prior  work  on  the  thermodynamic  properties  of  this  region  includes  the 
calculations  by  Mullins,   et  al .    (1962)  .     The  thermodynamic  properties  along 
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the  vapor  boundary  of  this  region  are  available  from  the  intersection  of  the 
32-terin  MBWR   (Weber,   1977b)   and  the  solid  vapor  pressure  curve  defining  the 
boundary  of  the  region.     The  elements  required  to  complete  the  description 
of  the  region  are  the  reference  point  for  all  properties,   the  vapor  pressure 
equation  of  the  solid,   the  density  of  the  solid,   and  the  heat  capacity  of 
the  solid. 

6.1    The  Reference  Point 

The  reference  point  for  this  region  is  defined  in  exactly  the  same  way 

as  for  hydrogen   (see  section  4.1).     The  heat  of  fusion  at  the  triple  point 

was  experimentally  measured  by  Giauque  and  Johnston   (1929)   to  be  106.3  cal/mol 

(444.8  J/mol) .     The  triple  point  temperature,   54.359  K,   and  the  melting 

pressure  equation  are  taken  from  Weber   (1977b).     The  volume  change  on  fusion 

3 

at  the  triple  point  becomes  approximately  0.94  0g  cm  /mol. 

6.2    Vapor  Pressure  of  the  Solid 

The  same  general  remarks  apply  to  oxygen  as  they  did  to  hydrogen.  The 
vapor  pressure  equation  is  of  the  Kirchhof f -Rankine  type 

log  P     =         1096^562485     _     2.025578307  log  T  +  28.35976524 

with    T    in  kelvin  and    F     in  mm  Hg.     The  constraints  used  are: 

1.  P^     =     0.0014451  atm,   T^     =     54.359  K, 

2.  |||       =     0.000482427  atm/K,  and 

3.  P     =     0.000019200  atm  at  44  K. 

Intercomparisons  with  the  available  experimental  vapor  pressures  of 
Aoyama  and  Kanda   (1936)   down  to  the  S-y  solid  transition  show  departures  of 
up  to  30%  in  pressure. 

6.3     Density  of  the  Solid 

Oxygen  is  an  interesting  solid  because  it  occurs  in  three  different 
solid  structures.     We  are  interested  here  only  in  the  so-called  y-solid 
which  exists  between  the  temperatures  of  43.801  K   (Kemp  and  Pickup,   1972)  and 
the  triple  point,   54.359  K.     For  this  range  of  temperatures  the  thermal 
expansion  measurements  of  Manzhelii,   et  al .    (1966)    the  absolute  measurement 
of  Tolkachev  and  Manzhelii   (1965)   and  the  x-ray  measurements  of  Barrett, 
et  al.    (1967)   are  in  reasonable  agreement,   i.e.,  within  0.6%.     We  represent 
the  experimental  results  with  a  parabolic  function  forced  to  pass  through 
the  designated  value  at  the  triple  point.     The  relation  is 
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=     V  =     23.2808187  -  0.06772142868  T  +  0.001339285715 


p     ,  .  ,  solid 
solid 

with  V     ,  .  -  in  cm'^/mol  and    T     in  K.     The  solid  densities  should  be  in  error 
soiia 

by  no  more  than  0.6%  in  the  range  of  temperatures  covered. 

6.4     Heat  Capacity  of  the  Solid 

The  experimental  measurements  on  the  y-solid  by  Giauque  and  Johnston  (1929) 
were  represented  by  Mullins,  et  al.    (1962)  with  a  simple  polynomial  which  is 
retained  here 

C^^„     =     16.908081  -  0.24181777  T  +  0.0024809089  T^ 
SAT 

with  Cg^rp  in  cal/mol   (4.184  J/cal)   and     T     in  kelvin.     Departures  of  the 
polynomial  from  the  measured  points  are  no  larger  than  the  stated  error  in 
the  measurements,   i.e.,   about  0.5%  in  Cg^^f 

6.5     Thermodynamic  Consistency  Checks 

These  checks  are  performed  as  described  in  section  4.5  for  hydrogen, 
except  that  here  the  temperature  limitation  is  that  of  the  B-y  transition, 
i.e.,   43.801  K.     Results  are  shown  in  table  3.     Departures  between  the 
different  ways  of  calculating  the  heat  of  sublimation  do  not  exceed  0.06%. 

Table  3.     Intercomparison  of  the  Keats  of  Sublimation  of  Oxygen 

AH     ,  ,  .      ^  .  AH     ,  ,  .      ,  . 

sublimation  sublimation 

T,   K  from  the  boundaries  from  the  C-C  Eq.  Difference 

J/mol  J/mol  J/mol 

44.000  8376.665  8376.109  .556 

45.000  8359.501  8359.234  .267 

46.000  8342.413       ,  8342.344  .069 

47.000  8325.380  8325.433  -.054 

48.000  8308.379                           -       8308.497  -.118 

49.000  8291.388  8291.527  -.139 

50.000  8274.383  3274.512  -.129 

51.000  8257.341  8257.443  -.102 

52.  000'  8240  .  238     .    "    •,  8240.305  -.  067 

53.000  8223.050  '  8223.083  -.033 

54  .  000     .        .      ■      8205.  751  8205.758  -.  007 

54.359  ■   '       8199.509  8199.509  defined  to  be  the 

same 
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7.      THE  COMPUTER  PROGRAMS  FOR  SLUSH  HYDROGEN  AND  SLUSH  OXYGEN 

The  phase  boundaries  and  representations  of  the  different  properties 
along  the  phase  boundaries  described  in  the  previous  sections  have  been 
combined  into  two  computer  programs.     The  programs  are  in  the  form  of 
subroutines : 

SLH2(TT,  QQ,  TP,  PSL,  DSL,  HSL,  SSL,  USL) 
SL02(TT,   QQ,   TP,   PSL,   DSL,   HSL,    SSL,  USL) 

where  the  variables  in  the  parameter  list  are  temperature,  quality,  triple 
point  key,  pressure,  density,  enthalpy,  entropy,  and  internal  energy. 

7 . 1  Input 

The  user  must  specify  temperature  and  quality.     He  may  specify  the 
triple  point  key  depending  on  the  temperature. 

Temperature .     The  input  temperature  in  kelvin  should  be  between  10  and 
22.60  K  for  hydrogen  and  between  44  and  58.2  K  for  oxygen.     The  programs 
will  return  answers  for  temperatures  less  than  the  lower  limit,  but  the 
values  returned  become  less  and  less  reliable  the  lower  the  temperature.  For 
temperatures  greater  than  the  upper  limit  an  error  message  is  printed:     "T  is 
too  high,  P  will  exceed  5000  psia." 

Quality.     A  value  for  quality  must  be  entered  such  that  0.0  <_  Q  <_  1.0. 
A  quality  of  0.0  will  cause  the  subroutine  to  return  values  on  the  denser 
side,   i.e.,   the  values  of  the  solid.     A  quality  of  1.0  causes  the  subroutine 
to  return  values  on  the  less  dense  side,   i.e.,  values  of  the  vapor  for  the 
solid-vapor  region   (temperatures  below  the  triple  point) ,  and  values  of  the 
liquid  for  the  solid-liquid  region   (temperatures  above  the  triple  point) .  If 
the  input  Q  is  fractional  the  subroutine  calculates  the  values  of  the 
properties  on  both  boundaries  of  the  region  and  then  performs  the  ratio 
according  to  the  equations  given  in  section  2.2. 

Triple  point  key.     For  the  special  case  when  the  input  temperature  is 
exactly  equal  to  the  triple  point  temperature,   13,800  K  for  hydrogen  and 
54.359  K  for  oxygen,   the  user  must  specify  through  the  triple  point  key,  TP, 
in  which  region  he  wants  the  subroutine  to  calculate  the  properties. 
TP  =  0.0  means  that  the  program  will  return  properties  from  the  solid-liquid 
region^  while  TP  =  1.0  insures  that  the  program  calculates  the  triple  point 
properties  in  the  solid-vapor  region. 

7.2  Output 

The  output  from  the  program  is  handled  through  the  parameter  list.  For 
a  given  input  of  T  and  Q    (TP  if  appropriate)    the  output  is  unique,  i.e., 
single-valued.     The  values  returned  are  the  slush  pressure  PSL  in  atm,  the 
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slush  density  DSL  in  mol/L,   the  slush  enthalpy  HSL  in  J/mol ,  the  slush 
entropy  SSL  in  J/mol-K,   and  the  slush  internal  energy  USL  in  J/mol. 

7.3     Flow  Charts 

The  slush  routines  are  detailed  in  the  flow  chart  of  figure  3.  The 
various  phase  boundaries  and  the  interiors  of  the  solid-liquid  and  solid-vapor 
regions  are  keyed  according  to  the  inset,  which  is  just  a  repeat  of  the  phase 
diagram  already  shown  in  figure  2. 

7.4     Program  Listings,  Hydrogen 

A  listing  of  the  hydrogen  program  is  given  in  Appendix  A.     The  slush 
hydrogen  routine  calls  on  the  program  deck  of  the  parent  32-term  MBWR  (Roder 
and  McCarty,   197  5) .     Those  routines  of  the  parent  deck  that  had  to  be 
modified  to  insure  consistency  around  the  triple  point  are  also  listed  in 
Appendix  A.     They  are: 

FUNCTION  PRESSM(TT) ,  and 
SUBROUTINE  RHOl (PP , DD , TT) . 

New  routines  developed  for  the  slush  application  and  listed  in  Appendix  A  are: 

FUNCTION  DPDTMP(TT),  and 

FUNCTION  FINDM(P,T,D) . 
The  first  is  simply  the  temperature  derivative  of  the  melting  pressure  curve, 
while  the  latter  is  an  additional  density  finder.     The  FINDM  routine  differs 
from  the  normal  FINDD  density  finder  in  that  it  requires  a  reasonable  trial 
density  in  its  calling  parameters. 

7,5     Program  Listings,  Oxygen 

A  listing  of  the  oxygen  program  is  given  in  Appendix  B.     The  slush 
oxygen  routine  calls  on  the  program  deck  of  the  parent  32-term  iMBWR  (Weber, 
1977b) .     Routines  not  available  with  the  parent  deck,  or  those  that  had  to 
be  modified  are  listed  in  Appendix  B  as  well.     They  are: 

FUNCTION  FINDT(P,D) 
SUBROUTINE  T1(PP,DD,TA) 
FUNCTION  PMELT(PP) 

FUNCTION  PRESSM(TT)  .V 
■   FUNCTION  DPDTMP(TT) 

7.6     Test  Programs,   Sample  Output 

A  sample  program  for  hydrogen  is  shown  in  Appendix  C ;   in  Appendix  D  for 
oxygen.     The  output  from  the  test  program  for  hydrogen  in  the  form  of  an 
abbreviated  table  of  values  with  particular  attention  to  the  triple  point  is 
given  in  table  4;   the  companion  output  for  oxygen  in  table  5. 
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Subroutine  SLH2  (TT,aQJP,PSL,OSL,HSL.SSL.USL) 
Subroutine  SLO2  (n.QQJP,PSL.DSL.HSL,SSL.USL) 


Enter  with  temperature  TT  and  quality  QQ 


QQ  =  0  Properties  on  the  denser  side 
QQ  =  1  Properties  on  the  less  dense  side 
QQ  Fractional  — ^  mixture  properties 


A  -  Solid  below  TP 
B  -  Vapor  below  TP 
C  -  Mixture  below  TP 


D  -  Liquid  above  TP 
E  -  Solid  above  TP 
F  -  Mixture  above  TP 


Special  Case,  Temp.  TT  =  Triple  Point 
TP=0.  for  Solid-Liquid  Properties 
TP=1.  for  Solid-Vapor  Properties 


^    Calculate  vapor  properties  below  T.P 


^    and  mixture  properties  if  appropriate 


Figure  3.     Flow  Chart  of  the  Subroutines  SLH2 , 


SLO, 
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The  sample  programs  are  intended  to  serve  as  test  decks  for  running  at 
other  installations.     They  include  a  call  to  every  possible  combination  of  T 
and  Q  in  the  regions  solid-liquid  and  solid-vapor.     To  round  out  the  list  of 
examples  the  liquid-vapor  region,  which  is  obtained  from  the  parent  32-term 
MBWR,  has  also  been  illustrated. 

The  astute  observer  will  note  that  the  various  properties  in  lines  4  and 
7  of  tables  4  and  5  and  also  lines  6  and  14  and  lines  9  and  13  should  be 
pairwise  identical,  except  of  course  for  Q  and  TP.     To  the  extent  that  these 
lines  are  not  identical  thermodynamic  inconsistencies  remain.     One  of  the 
inconsistencies  arises  from  the  fact  that  we  carry  10  digits  in  the 
coefficients  of  the  vapor  pressure  equations  and  the  melting  line.  The 
pressures  returned  from  each  of  the  different  functions  at  the  triple  point 
temperature  are  very  nearly  the  same,  however,  they  are  not  identical  in 
every  bit  of  the  computer  word.     Entering  the  32-term  equation  of  state  with 
slightly  different  values  of  the  pressure  results  in  slightly  different 
values  of  vapor  densities  returned.     The  disagreement  in  entropy  between  lines 
6  and  14  arises  from  an  inconsistency  inherent  in  the  equation  of  state.  This 
is  shown  simply  by  the  fact  that  the  AH  between  lines  13  and  14  of  tables  4 
and  5  divided  by  the  triple  point  temperature  is  not  equal  to  the  AS  betv/een 
lines  13  and  14.     Several  other  discrepancies  are  attributable  to  rounding 
errors.  ;  ■ 

7.7     Unit  Conversions 

Common  conversion  factors  for  use  with  hydrogen  are  given  in  Appendix  E, 
while  those  for  use  with  oxygen  are  shown  in  Appendix  F. 

8.      ERROR  ESTIMATES,   FURTHER  WORK 

Whenever  possible  estimates  of  uncertainty  were  given  in  the  sections 
that  described  the  fitting  or  representation  of  a  particular  piece  of 
experimental  data.     We  summarize  these  uncertainties  as  maximum  errors  in 
the  different  properties  used  in  the  computer  programs  in  table  5. 
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Table  6.     Estimates  of  Maximum  Uncertainties 


Variable 

temperature 

pressure 

density 


Hydrogen 


Oxygen 


0.01  K 


0.01  K 


1%  at  the  triple  point  1%  at  the  triple  point 

increasing  to  200%  at  very         increasing  to  31%  at  44  K. 
low  temperatures.     Negligible  Negligible  for  pressures  above 


for  pressures  above  1  atm. 

0.2%  at  temperatures  above 
T.P.     0.4%  at  the  triple 
point  increasing  to  1.5%  at 
4  K. 


1  atm. 

0.2%  at  temperatures  above 
T.P.     0.6%  for  temperatures 
below  the  triple  point  and 
above  the  6-Y  transition. 


enthalpy, 

internal 
energy 


entropy 


1.2  J/mol  for  solid-liquid 


1.5  J/mol  for  solid-vapor 


0.08  J/mol-K  for  solid- 
liquid 

.146  J/mol-K  for  solid-vapor 


3  J/mol  at  the  triple  point 
increasing  to  27  J/mol  at  the 
highest  pressures. 

3  J/mol  at  the  triple  point 
increasing  to  6  J/mol  at  the 
lowest  temperature. 

0.05  J/mol-K  at  the  triple  point 
increasing  to  0.46  J/mol-K  at 
the  highest  pressures. 

0.05  J/mol-K  at  the  triple  point 
increasing  to  0.13  J/mol-K  at 
the  lowest  temperature. 


Entropy  errors  are  approximated  by  considering  enthalpy  errors  divided 
by  the  temperature  and  vice  versa.     For  hydrogen  the  entropy  departure  for 
the  solid-vapor  region  is  taken  from  the  closed  loop  third  law  calculation 
of  Roder,  et  al.    (1965) .     For  the  solid-liquid  region  we  use  the  experimental 
estimate  of  1%  in  measured  heats  of  fusion  published  by  Dwyer,  et  al.    (196  5b) . 
For  oxygen  we  use  a  value  of  3  J/mol  as  error  estimate  in  enthalpy  (Weber, 
1977b)   and  add  to  this  either  6%  of  the  heat  of  fusion  to  get  the  solid-liquid 
departures,  or  the  enthalpy  equivalent  of  0.5%  of  the  values  of  C^^rp  to  get 
the  solid-vapor  deviations.     The  rather  large  discrepancies  described  in  an 
earlier  version  of  this  report   (letter  report,  June  1977)   have  indeed  been 
resolved,   i.e.,   reduced  considerably,  by  using  a  PVT  surface  which  is 
consistent  with  the  melting  line  and  the  triple  point. 

The  work  described  here  has  led  indirectly  to  a  review  of  the  molar 
volume  of  solid  oxygen   (Roder,   1977).     One  conclusion  of  that  review  is  that 
the  volume  change  on  fusion  at  the  triple  point  is  uncertain  by  about  0.6%. 
Further  experimental  work  to  accurately  measure  this  volume  change  at 
essentially  zero  pressure  is  recommended.     Additional  experimental  work  on 
oxygen  may  be  required  particularly  if  large  scale  use  of  slush  oxygen  and 
hydrogen  is  likely.     The  sponsor  should  arrange  to  have  the  dielectric 
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constant  of  oxygen  reviewed  and  measured  for  the  solid  since  this  property 
will  be  required  for  quantity  gaging.     In  section  5  we  indicate  that  the 
thermodynamic  properties  for  the  solid-liquid  region  for  oxygen  had  to  be 
estimated.     Some  experimental  measurements  are  required  to  define  the 
properties  adequately.     We  suggest  that  measurements  of  the  heat  of  fusion 
or  the  volume  change  on  fusion  be  conducted  for  oxygen  at  elevated  pressures. 
An  alternative  possibility  is  to  measure  heat  capacity  values  for  oxygen  in 
the  solid  state.     Finally,  a  review  of  Prandtl  and  Grashoff  numbers  should 
be  made,  to  include  values  of  C  ,  C     thermal  conductivity  and  viscosity. 

P  V 

Several  of  these  properties  are  presently  not  defined  for  a  slush  regime,  yet 
they  will  be  needed  to  define  heat  transfer  characteristics. 

9.  SUMMARY 

Computer  programs  to  describe  the  thermodynamic  properties  of  slush 

hydrogen  and  slush  oxygen  have  been  developed.     The  range  of  temperature  for 

these  programs  is  from  several  degrees  below  the  triple  point  to  several 

degrees  above  the  triple  point.     The  upper  temperature  limit  coincides  with 

2 

a  pressure  limit  of  about  5000  psia   (34  MN/m  ) .     In  the  case  of  oxygen  the 
experimental  data  are  very  limited  and  the  properties  calculated  are  in  the 
nature  of  a  "best  estimate."     The  equilibrium  and  thermodynamic  properties 
included  are  temperature,   pressure,   density,   enthalpy,   entropy,  internal 
energy,   and  quality.     Input  to  the  programs  is  temperature  and  quality,  the 
remaining  variables  are  outputs. 
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APPENDIX  A.      PROGRAM  LISTINGS  HYDROGEN 


SU3Rl)UTINE  SLH2(TT,UQ,TP,PSL,0SL,HSL,53L,JSL) 

THERMOO YMmMIC  PROPERTIES  FOR  SOLIO-VAPOR  OR  SOLID-LIQUID  REGIONS 
TEMPERATURE   TT    MUST    BE   BELOW   22.6  K 

QUALITY   QQ,    0.    IMPLIES   PROPERTIES   ON   THE  OESSER  SIDE 

1.    IMPLIES   PROPERTIES   CN  THE  LESS   DENSE  SIDE 
FkACTIONAL   QQ  BETWEEN   J.    AND   1.  ?ERMISSA3LE 
OTHER   VALUES    OF  QQ  INVALID 

FOR   THE   SPECIAL  CASE   OF   T=TRIPL£   POINT   TP=0.    MEANS  SOLID-LIQUID 

TP=1.    MEANS  SOLID-VAPOR 

DATA    (VPl=-90  .775  66  9^9)f{VP2  =  2.^69<330  9'+u),IIVP3='*.00  9a5735i4) 
DATA    tCPl  =  -i3l.)»(CP2  =  d7.75),(CP3  =  -26.95833333),  (CPi*=k,25)  , 
1    (CP5=-0. 0Hlo66b6c67) 
DATA    (Vl  =  22.7  79  63'+74)  ,(V2  =  -O.O3l),(V3  =  J.O05) 

DATA    (TTP  =  13. 800),  tTUP  =  22. 63),  (H  TP  =  -7 '♦3.206291)  ,  (  ST  P  =  l  .49995  90  ) 
T=TT 

Q=aQ 

PSL  =  OSL=iHSL=SSL=JSL=  u  . 
IF ( T . LT .TUP  )    GO   TO  1 

PRINT  2 

FORMAT ( 5 X , 3dHT    TOO   HIGH,    P   WILL    EXCEED   5000  PSIA) 
RETURN 

IF(u.GE.U.    .AND.   U.LE.l.j)    GO   TO  5 

PRINT  4 

FORMAT ( 5 X , ItHG   OUT   OF  RANGE) 
RETURN 

IF(T-TTP)  6,7,6 

ENTER  HERE   FOR  TRIPLE  POINT 
IF ( TP.EQ. J . j)    GO   TO  a 

ENTER  HERE   FOR   SOLIu-\^APOR  PROPERTIES 
P  SL  =  E  XP  (  VP  1/T  +  VP2*AL0G  (T  )  +\/P3) /76a  . 
0P0T  =  PSL»C-VPi/T»*2«- VP2/T) 
VSCL= V1*V2*T+V3*T**2 
OSL=130G./VSOL 

HSL  =  HTP-(CP1*(TTP-T)4-CP2/2.*(TTP»*2-T*»2)+C  =  3/3.*(TTP»*3-T»*3) 
1    <-CP4/4.*(:TTP»»'*-T»»u)  +  GP5/5.*(TTP»»5-T»*5)  )  /  lOOO. 

SSL  =  ST3-(CP1*{AL0&(TTP)-ALGG(T))+CP2*(TTP-T)4-CP3/2.*(TTP**2-T**2) 
1   <-GP4/3.*(TTP»*3-T»*3)+CP5/4.»(TTP**4-T*»i^))/1000  . 

USL  =  H5L-PSL/LJSL»i:i.  325 
USOL=USL 

IF (Q. EQ  .  G  .  L  )    RETURN  ,  - 

CONTINUE   hcRt.   FOR   VAPOR  PROPERTIES 
OVAP=FINDM(FSL, T, C. u6) 
VVaP=103u ./DVAP 

DEL V=VVAP-VSOL  .  .  " 

0£LM=CELV*T»GFDT* G. 1G1325 

GELS=OELH/T 

UVAP=HSL+DELH-PSL/(1GGO./(OELV+1033./DSL))»1Q1.325 

DELU=UVAP-USOL 

IF (Q. NE . 1 .  ;  )    GO  TO  9 

GSL  =  lQ03./(LiELV  +  li3GC./C5L) 

HSL  =HSL  +u£Lrl 

UGL  =  HSL-PSL/GSL»-iGl.325 

SSL=5SL+DELS 

RETURN 

CONTINUE    HERE   WiTh    FRACTIONAL  Q 

0  3L  =  1'jGG./(Q*l)ELV  +  1G:K./DSL) 
HSL  =  u*OELH-t-HSL 
SSL  =  Q*QclLS+SSL 

2  6 


RETURN 

ENTER  -^LRl  FCR  SOLlD-LliUlJ  ^RO^ERTIES 
PSL=?RE3SM ( T ) 

dfgt=up:tmf (T J 

3SL=FINDM  (F3L  ,  T  ,    3  .  ) 
H5L=ENThAL(PSL,3SL,T1 
USL=HSL--SL/aSL*l01. 325 
JLIG=USL 

SSL=ENTkCP(3SL,T) 
IF(G.E3.1.a)  RETUkN 

CONTINUE   rMRE   FCR  SOLID  PROPERTIES 
JELH=(25.0^o9  3l5'+♦C.0'♦i♦lH9♦  =  SL)♦'^•.18i^ 
DELV  =  OE^..-t/T/OPGT/0.a2!»2i79/H.16-^ 
CELS=OELH/T 

OSL  =  ii:uu./{ljOO./CSL-G£LV) 

USCL=iSL-OELH-PSL/DSL*lui. 325 

DELU=ULIG-u3Cl 

IFCQ.NE.b..)    3  0   TO  13 

HSL  =  HSi.  -GELH 

S5l=SSL-CEl3 

USL=USOL 

RETURN 

CONTINUE   t^lRE.    rilJ-i    FRACTIONAL  G 

OSL  =  130:./(.**DELV  +  1CJC./:3L) 

rtSL=rtSL-G*OELH 

SSL=S5L-G*:EL5 

USL=USL-G»GELU 

RETURN 

END 


FUNCTION  FR£SSM(TT) 

CALCULATES   '-.ElTING   PRESSURE   FRO'I         InPjT  TEMPERATURE 
T=TT-t-j,JC3 

PS=.Gc95+(T-13.3:3)*3:.33i2»EX=  (-5.d93/T)+(T-13.3Q3)*2.0»T/3.0 

PRE5S.M=?S 

RETURN 

END 


SUSRO'JTINE  RHul(PP,GD,TT) 
1ST   CUT    AT    RHO   FROM    F  =  A-t-3*T 

REALLt    an    iteration,    BUT    IT   MAY    EE   SMA^L   ENOUGH    AND    FdST  ENOUGH 
P    IN    ATM,    T    IN   <,    RHO   IN  MGLES/lITER 
DIMENSION   ncH0(i*3)  ,A  (i*3)  ,5  (-3) 

DATA((RHC(I),I=l,'.3)=l.,2.,3.,u.,5.,b.,7.,3.,9.,iO.,ll.,l2.,13.,lit 
1., 15. ,16. ,17. ,13. ,19., 20. ,21.  ,22. ,23.  ,2*.. ,25. ,2b. ,27. ,23. ,29.  ,3  3., 
25l.,32.,33.,3'*.,35.,3c.,37.,3c,,39.,'*0.,'+l.,'»2.,'+3.) 

DAT4((A(I),I=l,^3)=-.Z9599325  7i*l,-i.l335'+55l75,-2.i*35l37S102,-^.3: 
ItrlCSg  31  7, -..;  =  9t5  5393, -i.2559135Hj5, -12. 33c2b6955, -15. 733961531,-1 
29. 57337;  659, -23. o3i:-*  1295, -23.036^58235, -32. 77255895  1  ,-37. 72663671 
39, -•*2.:;«»£371352, -2--325>*l,-5'f.i9i, '+6  51^2, -60. 237558136, -66. 7267 
•♦'3 >»77l, -73. 516513375, -  .■)i.6-f3765-*19, -33.  j7990  =  s9'*, -95. 751501525,-133 
5  . 6  3  96  :;  ^1,-1 11.  72  777  3*.  2  , -119.  799  7- 17  2  , -127.  786  579  , -  135.  9^-257  76  1  ,-1 
6-»3.3333i332, -151. 265■+9•:)D-^, -157.  996  3  8i+63  , -163. 7133357, -167. 27873589 
7,-169.7l92<:3;3,-169.769939l-n-lo6.3;7233-*5,-16w.22  9  7o-f5'*,-l-*9.6;3l 
3  3c3,-l3-.=;3w;21o,-116.1l756-f;i,-9-t.293ll7563,-67,52-f2-*7269,-3c.7* 
923953  2-*, -1.  32  9557^977) 

DATA((3(I),I=l,-.3)=.  3  35j  5  92335-9,.  1^  575l-:3596,.2717Eo7:*97  9,.3  72  9  76 
113  72:*,.ii.7952911-kc5,.5^:14l2766l5,.7j36-»366952,.33l2'*6-i>+lH2,.95923U-f 
291  l-,l,,923iu2o55,1.23223677il, 1.3776-55569,  1  .52  9'*  993193, 1. 6533212 
3266, 1.35-32-(-^9l3, 2.  33;  335-23, G.2l5i;C6l  3  9, 2. -'ll2o7-*7  3  5, 2. 6l95u3366 
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45  ,Z  .3^0  7G't32l2,  3.  J756'f5J8  76,3.  3  250  4a55'+  ,3.  5895  52500  8,  3  .de968g35i+l, 
5^+.  165  85  93639, i*. '+73296  l9'ti+,<+. 820^02  20  69,  5.  181*0092761  ,  5  .  5696'+6053  ,5  . 
697725  2^+ 362,  6.  it057l7553l, 6.  321 '♦2  ^22  39,  7.  2730  6  97  825,  7  , 7 1+ 046 1 65  3  U  ,  3  .  1 
7  9  725755a2,d.6358'46a'+65,9.J'+73286  C7&,9.-»26  0953807,  9.  81 32535  0^+5  , 10  .  2 
822281 36 3,  1J.03Q12934, 11.  05913823  , 11. 56530719'*) 
10  P=PP 

T  =  TT 

IP  =  1 

PLO=0. 

C  CHECK   INPUT   DATA   FOR  RANGE 

C         INPUT   RESTRICTIONS  REMOVED 
C  IF{P.GT.35G.>    GO  TO  15 

C  IFiT.GT. 300.1)    GO   TO  15 

C  GO   TO  13 

C     15   PRINT  16 

G     16   F0RMAT(2'+H   INPUT   GATA   OUT   OF  RANGE) 
C  0  0=u.Q 

C  RETURN 

PTP=\/PN(l3.a) 

IF ( P. GT . PTP )    GO  TO   1 8 

IFCT.GT .13.3)    GO   TO  16 

oc=Q.caaoi 

RETURN  - 
C  CHECK   MELTING  LINE 

13    IF(  T.GT  , (+3  .  60  0)    GO   TO   23  . 

PS=='RESSM  (  T  ) 

IF(P-PS)  23,21,19 
19   PRINT  22 

22  F0RMAT(26m    INPUT   CONDITIONS   IN  3CLI0) 
21  30=39.3 

RETURN 

G  PhASEFINOER  •  . 

23  IF( T.LT . 32. 938)    GO   TO  27 
PFHASE=-5G . 6u 02 +1 .92  0 888*T 
IF (P.GT.PPriASE)  IP=15 

GO   TO  36 
27    P\/AP  =  VFN(T) 

IF(P-F\/AP)  36,  32,35 

32  PRINT  33 

33  FORMAT(53H    INPUT   PLACES    YOU   EXACTLY   ON   THE    VAPOR   PRESSURE  CURVE) 
3J=0SATL  (T) 

RETURN  .     ■  . 

3-:   OEfJL  =  DSATL  (  T) 
IF  =  QErJL 

C  START   TA3LE   LOOKUP  HERE 

3d    CONT INUE 

00   <+0    I  =  IP,h3 

PCALC  =  A  (  I)  +d( I) *T 

IF(P.LT.PCALC)    GO   TO  kl 

PLO=PCALC 
"^J    CON  TI  NUE 

IF  (  P  .  LT  .  PCALC)    I  =  -»3 

IF(P.LT.FCmlG)    go   TO  i.1 
C  -PRINT  ^7 

C     -7    FGRMAT(35h   HIGH    DENSITY,    OUT   OF    RANGE   FOR  NOW) 

0  L  =  H  •+  , 

RETURN 
Hi  CCNTIfJUE 

00  5Hj=i,ij 
RHOF=J 

RHOF  =  ,-^H0F/ IC  . 

IF(I.EQ.l)    GO   TO  5u 

AA  =  A  (I-i 3  +  (A  I  I ) -A  ( I- 1) )  tRH OF 
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33=3(I-i) *(3(I)-B(I-1) )*RHGF 

GO  TO  52 
33  AA=A(I)*RHGF 

33=8(I)»kH0F 
52  PCALC=AA-t-3o*T 

IF(P.LT.PCALC)    GO   TO  55 

54  ?LO=PGALC 

55  F'^AC=  (P-PLu  )  /  (PCALC-PLO) /lij . 
00  =1-1 

uu     =00      +KhOF-0. 1+FRAC 

RtTURN 

ENTRY  Ti 

FIRST  GUESS  FOR  TEMPERATURE  ITERATION  OF  FINOT 
P  =  PP 
0  =  00 

00  60  1=1.^3 

IF(0.LT.RH0(I) )    GO   TO  61 

c3  CONTINUE 
I='*3 

61  FPAC=0-RHO ( I-l) 
IF(I.EQ.i)    GO   TO  63 
AA=a(I-l)<-(A(I)-A(I-l))  ♦FRAG 
3B=B(I-1J*-  (3{I)-a(I-l))*FRAC 

GO  TO  c2 
63  FRAC=0 

AA=A(I)»FRAC 
33=3  (  I)  ♦FRAG 

62  TT=(P-AA)/33 
END 


FUNCTION  OPaTMP(TT) 

NEW   FOR   HYORGGEN*    .BELTING   PRESSURE  DERlVAriy/E 
OATA  (TTP=l3.8u3),(Al=3Q.33l2),{AL=5.b^3) 
T  =  TT  +  C. GO  3 
A2=2./3. 

EXPP=Ai»EXP (-AL/T ) 

0PQTMF=(T-TTP)*(EXPP*AL/T**2+A2)+EXFP+A2*T 

RETURN 

ENO 


FUNCTIUN   FINO  .'^(P,T,0) 
C  CHANGEO,    TRIAL   OENSITY    VIA   SUBROUTINE    RHOl/Tl,    Zk   FES  1975 

Z  CHANGED    TO    Aui-lIT    A   TRIAL   DENSITY    FOR    SLUSH    USE      15    DEC  76 

COM VON/D«TA/G,R, GAMMA, VPfuTP 

DIMENSION   G( 32)  ,Vr ( 3) 

TYPE    DOUoLEGi R, GAMMA 

TT  =  T 

IF ( a. NE  .  u  .  )    GO   TO  2 
CALL  RhG1(P,j,TT) 

:  DO=D 

DO   IJ   1  =  1, 5j 

C^LL   PRESS ( PP , DO, IT) 

.P  2  =  P  P 

IF(A35    {  P-P2 )  -1  .E-7*c ) 20  ,  20  , 1 
i    C-LL  DPDD(PP,OD,TT) 

DP  =  PP 

CjRR=  (P2-P) /L.t- 

J  =  uC 

IF(A3S    (CORR) -1 .E-7» D ) 20 , 2u , lu 
1^  QO=Ju-CGRR 
2  :    F  I    D   M=  -]  J 

ENO 
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APPENDIX  B.      PROGRAM  LISTINGS  OXYGEN 


SUBROUTINE   SLO 2 ( TT , QQ , TP , PSL , OSL ♦ HSL t S SL , US L ) 

C  THERMODYMAMIC  PROPERTIES  FOR  SOLID-VAPOR  OR  SOLIO-LIQUIO  REGIONS 
C  ■        TEMPERATURE   TT   MUST   9E   BELOW  58.2  K 

C  QUALITY   QQ,    0.    IMPLIES   PROPERTIES   ON   THE  DENSER  SIDE 
C  1.    IMPLIES   PROPERTIES   ON   THE  LESS   DENSE  SIDE 

C  FRACTIONAL  QQ   BETWEEN   0.    AND   1.  PERMISSA3LE 

C  OTHER   VALUES   OF  QQ  INVALID 

C  FOR   THE   SPECIAL   CASE   OF   T=TRIPLE   POINT   TP=0.    MEANS  SOLID-LIQUID 
C  TP=1.    MEANS  SOLIO-VAPOR 

DATA   (VPl  =  -1096.562i»85)  ,  (VP2=-2.  025578307)  ,(VP3  =  28.3597652'*) 
DATA    (  CP1  =  16.  908081)  ,  (CP2  =  -0.  2 '♦18 1777)  ,  ( CP3=  C  .  0  02i»8  09  0  89) 
DATA    (Vl  =  23. 2803137)  ,  (V2  =  -0  . 067721'+2868)  ,  (V3  =  0. 001339285715) 
DATA    (TTP  =  5't.  359)  ,  (TUP=5  8.20)  ,  ( H  TP  =  -66  25  ,  9  7  72  6  )  ,  (  STP  =  5  8  .  928i*7  8 ) 
T  =  TT 
Q  =  QQ 

PSL=OSL=HSL=SSL=USL=0 . 
IF (T.LT . TUP)    GO   TO  1 
PRINT  2 

2   FORMAT (5X,35HT   TOO  HIGH,   P   WILL   EXCEED   5000  PSIA) 
RETURN 

1   IF(Q,GE.O.    .AND.   Q.LE.i.O)    GO  TO  5 

PRINT  k 

U   F0RMAT(5X,14HQ   OUT   OF  RANGE) 
RETURN 

5  IF(T-TTP)  6,7,8 

C  ENTER  HERE   FOR  TRIPLE  POINT 

7   IF(TP.EQ.O.O)    GO  TO  8 
C  ENTER  HERE   FOR   SOLIO-VAFCR  PROPERTIES 

6  PSL=EXP ( VP1/T+VP2*AL0G (T) +VP3) /760 . 
DPOT=PSL» J-VP1/T»*2+VP2/T ) 

VSCL  =  V1  +  V2*T4-V3*T»*2 
DSL=1COO./VSOL 

HSL=HTP-(CPl*TTP+CP2*TTP*»2/2.+CP3*TTP*»3/3.-(CPl*T+CP2»T»*2/2. 
1   +CP3»T»»3/3.  )  )  18^ 

SSL  =  STP- (CP1*AL0G (TTP)  +CP 2 *TT P+CP 3 *TT P ** 2/2 . - ( CPl *  A LOG ( T) ♦CP2*T  + 
1   CP3*T»*2/2.)  )*it.l8i+ 

USL=HSL-PSL/OSL'101 . 325 

USOL=USL 

IF(Q.EQ.G.G)  RETURN 
C  ,:     ,   CONTINUE   HERE   FOR   VAPOR  PROPERTIES 
OVAP=FINDM« PSL,T ,0. 00  0  6) 
VVAP=1C00 ./OVAP 
OELV  =  V  VAP-VSOL 
DELH=OEL V*T* DP DT»0. 101325 
OELS=OELH/T 

UVA?=HSL+DELH-PSL/ (IQOC./ (DEL V^lOQO. /DSD) •101.325 
DELU=UVAP-USOL 
IF (Q.NE . 1. 0 )    GO   TO  9 
■OSL=1QCO  . /(OELVf  lO'lQ  ./OSL) 
HSL=HSL+OELH 
US L= HSL -PSL /OSL* 101. 32  5 
SSL=SSL+OELS 
RETURN 

C  CONTINUE    HERE    WITH    FRACTIONAL  Q 

9  OSL=10C0./(Q»OELV+1000./OSL) 
HSL=Q»0ELH4-HSL 
SSL=G*CELS+5SL  . 
•   USL  =  Q*0£LU4-USL 


30 


RETURN 

C  ENTER   HERE  FOR  SOLIO-LIQUIO  PROPERTIES 

8  PSL=PRESSM( T) 
0POT=0POTMP (T) 
OSL  =  FINOM(PSL,  T,i»0  .87) 
HSL=ENTHALf PSL ,OSL,T) 
USL  =  HSL-PSL/OSL*101.  32  5 
ULIQ=USL 

SSL=ENTROPCOSL,T) 
IFfQ.EQ.l.O)  RETURN 
C  CONTINUE  HERE  FOR  SOLID  PROPERTIES 

OELH=(  10  6.3+0.01126*PSL)*«+.13'+ 
OEL\/  =  OELH/T/0PDT/0.  0  2i»217g/i»,18t* 
DELS=OELH/T 

OSL=10  0  0./(10  0  0./DSL-OELV) 

USOL=HSL-OELH-PSL/DSL»101.325 

OELU=ULIQ-USOL 

IFro.NE.O.O)   GO  TO  10 

HSL=HSL-OELH 

SSL=SSL-DELS 

USL=USOL 

RETURN 

C  CONTINUE  HERE  WITH  FRACTIONAL  Q 

10  OSL=  100  0  ./(Q*DEL\/  +  100  0  ./DSL) 
HSL  =  >-SL-Q*OELH 
SSL=SSL-Q»DELS 
USL=USL-Q*OELU 
RETURN 
ENO 


FUNCTION   FIND  T{P»0) 
C  NEW  FEB  1975 

C  COPIED  FROM  HYDROGEN   FOR  USE  WITH  OXYGEN 

COHMON/OATA/G,R, GAMMA, VP » DTP 
DIMENSION  G(32)»VP(8) 
TYPE   DOUBLEGtR, GAMMA 
pp  =  p 

D0=0 

C  USES   A   FIRST   GUESS   IN   TEMPERATURE   FROM  Tl 

CALL  T1(PP,0D,TA) 
TT  =  TA 

DO   10  1=1,10 

CALL  PRESS(PP,OD,TT) 

P2  =  Po 

IF(A8S    (P-P2) -l.E-7*P) 20,20,1 
1   CALL   DPOT (PP.DO, TT) 
DP  =  PP 

CORR= {P2-P) /DP 

IF(A8S    (CORR3-1.E-5  520,20,10 
13  TT=TT-CORR 
20   FIND  T=TT 

RETURN 

•END 


SUBROUTINE  Tl(Po,DC,TA) 
C  SPECIAL    SETUP   FO"   TANK-P  OXYGEN 

TA=99. 

RETURN 
ENO 
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FUNCTION  PMELKPP) 

INVERSE  MELTING  PRESSURE  FOR  02,  INPUT  P,  ATM  OUTPUT  T,K 
NEW  WEBER  REFERENCES 

DATA    (P0  =  26  37.33),(CC  =  1.76  9)t<TTP  =  5i+.359),  ( P  TP  =  0  .  0  0  l^+^S  1) 

BRA=(PP-PTP  )/PO+1.0 

TEMP=ALOG«BRA) /CC 

PMELT=EXP (TEMP) *TTP 

RETURN 

END 


FUNCTION  PRESSM(TT) 

MELTING   PRESSURE  FOR   02,    INPUT   T,K  OUTPUT  P,ATM 
NEW  WEBER  REFERENCES 

OAT 4   (P0  =  2637.33)  ,  (CC  =  1.76  9)  ,  { TT PsS'* .  3 59 )  ,(PTP  =  0.00l'ti*51) 

PRESSM=PTP*PO» ( (TT/TTP) **CC-1.0) 

RETURN 

END 


FUNCTION  DPDTMP(TT) 

DERIVATIVE  OF   MELTING   CURVE  FOR  OXYGEN 
NEW   WEBER  REFERENCES 

DATA    (PO  =  26  37.33),(CC  =  1.769),(TTP  =  5^.359),(PTP  =  0.001'+it51) 

OPQTMP=PO»CC/TTP*»CC*TT** (CC-1.0) 

RETURN 

END 


APPENDIX  C.      TEST  PROGRAM  HYDROGEN 


PkCGRAM    values ( input , OUTPUT) 
C  THIS   PROGRAM   TESTS   THE    SLUSH   Hr3F.JG£N  SUBROUTINE 

CALL  OATAPHZ 

C  C0M3lrjATIG,NS   OF   SOLID    -  VAPOR 

PRINT  1 
T  =  12. 

a=G . 

CALL    SLH2 (T ,u» 0 . f PSL ,OSL , HSLf SSL  ,USL) 
PRINT   2,T,Q,        PSL , DSL ,HSL ,SSL, USL 
Q  =  J  .  5 

CALL  SLH2(T,Q,a.»PSL,CSL,HSL,SSL,USL} 
PRINT   3,T,U,        PSLt DSL ,H3L ,SSL, USL 
Q  =  l. 

CALL  SLH2(!T,a,a.tFSL»0SL,HSL,SSL,USL) 
PRINT  •♦,TtCi»        PSL, DSL, HSL, SSL, USL 
T=13. doa 
TF  =  1. 

a=!j . 

CALL   5LH2(T,a,TP,PSL, DSL,HSL,SSL ,USL) 
PRINT   5, T,vi, TP, PSL, USL,  HSL, SSL, USL 
Q  =  :.5 

CALL   SL H2 (T , a, TP, PSL , DSL , HSL, SSL  ,USL} 
PRINT   b,T,U, TP, PSL, DSL, HSL, SSL, USL 
Q=l. 

CALL  SLH2(T,a,TP,PSL,CSL,HSL,SSL,USL) 
PRINT   7, T,u, TP, PSL, DSL, HSL, SSL, USL 

C  CUMBlrjATIONS  OF   SOLID   -  LIQUID 

PRINT  i 
T  =  1  3  .  c  G  J 
TF  =  u. 
Q=G  . 

CALL    -,Ln2  {  T  ,  J  ,  TP,  FSL  ,  CSL  ,  HSL  ,  SSu  ,USLJ 
PRINT  :i,T,Q,TP,PSL,DSL,HSL,SSL,USL 
Q  =  J  .5 

CALL  3LH2(T,U,TP,FSL,CSl,H3L,SSL,USL) 
PRINT1u,T,u,TF,PSL,DSL,HSl,SSL,USL 

a=i. 

CALL    SL  H2  CT  ,  a,  TP,  ?:.L  ,  DS  L,  HSL,  SSL  ,  USL) 
PKlNTli,T,0,TP,FSL,CSL,HSL,5SL,USL 
T  =  15. 

a=c. 

CALL  iLHZtT,a,u.,l-SL,DSL,HSL,S3L,U3L!; 
'JRINTi2,T,U,        PSL, CSL, HSL, SSL, USL 
i  =  J  .  5 

CALL    SlH2  CT,Q,j  ,  ,FSi_  ,  CSL,  -'SLiSSl  ,USL3 
ORI NT  ii  ,  T  ,  :  ,        ?S  L , D  SL , HS  L  ,  SSL , U  SL 
J  =  1  • 

CALL   SLH2!;T,G,j.  ,^SL,GSL,HSl,33L,U3L) 
P^I\T1'*,T,Q,        PSL, DSL, H5u, SSL,  USL 

Z  CC.13I;.- TI  ji<5   oF   LI'iUIQ   -    ;APOk   FROi-^   T^E    32  TER^. 

PRINT  15 
T  =  l  3.  -i:  j 

?SL  =  ;r-j(T) 
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3SL  =  FINl)M(PSL, T, Jtt,  ) 
HSL=£NTHAL(PSL,3SL,T) 
USL=HSL-F5L/0SL*1Q1. 32  5 
SSL  =  ENTROP (OSL  ,  T) 

PRlNTlb,T,G,  PSL,uSL,H5L,S5L,USL 
Q=l. 

a5L  =  FINLM(PSL,T,0.a6) 

HSL=ENThAL (PSL » OSL , T )  n 

USL=HSL-PSL/DSL»lti. 325 

SSL  =  eNTROP  (  QSL  ,  T) 

PRrNTl7,T,C,  FSL,u5L,HSL,3SL,USL 

T=i5. aa J 

Q  =  j. 

PSL  =  \/PN{T) 

L)SL=FlNLiM<FSL,Tt38.  ) 
HSL=ENTHAL(PSL,OSL,T) 
USL=H5L-P5L/0SL*lLi.  3  25 
SSL=ENTROP(DSL,T) 

PRINTlbtTtu,  PSL,u3L,HSL,SSL,USL 
Q=l. 

OSL=Fir.!UM  (PSL  ,1,0.^6) 
HSL  =  ENTHAL (PSLi OSL , T  ) 
U5L=HSL-PSL/QiL»i01. 325 
SSL=£NTROP ( J5L , T) 
PRINTi7,T,U»        PSL, QSL , HSl  ,SSL 


7 

■3 
i  J 
11 
12 
13 


FOR.IAT 
1  ♦ 

FORMAT 
1    F  1  3 .  o 

FORMAT 
1    Fl 3. S 

FORMAT 
1    Fl 3. 3 

FORMAT 
1    Fl 3. d 

FORMAT 
i    F13. 3 

FORMAT 
1  F13.d 

FORMAT 
1  * 

FORMAT 
1   F13.  8 

FORMAT 
1  i 

FORMAT 
1  F13.6 

FORMAT 
1  F13.2 

FORMAT 
1    F  1  3  .  d 

FOR.IAT 
1  F  1  3.  0 
'FORMA  T 
1  ♦ 

FORMAT 
1    Fl 3.  3 

FORMAT 
1    F 1 3.  9 

E  NL 


///,52H1SGLIQ-VAP0R 
P,    ATM  RmO,  MOL/i 


Ig  X 
FlH 

iux 

Fl^i 
1  J  X 
Fl-f 
1  u  X 
Fli+ 
lOX 
Fl-^ 
1  u  X 
F  1^+ 
/// 
^  » 
1  U  X 
F11+ 
1  b  X 
Fl^ 
1  C  X 

'-"iH 

1  G  X 

Fm 
i  u  X 
Fl^ 
1  0  X 
FlH 
//  / 
?  , 
1  j  X 
Flw 
lb  X 
Fl^ 


F3.3,2<+H  SOLID, 

b,  Fll. 3 , F IC . 3 , Fll , 3 
Fo.3,2-fH    MIX,    J  IS 
6, Fll. 3,F1C.3,F11.3 
F3.3,2itH  VAPOR, 

6,  Fll.  3, F IC . 3  ,  Fll .  3 
Fd.3,2--H  SOLID, 
6, Fll. 3, FlC. 3,  Fll, 3 
Fj.i,2HH    MIX,    Q  IS 
o,Fll.3,Fli:i.3,  Fll .  3 
F6.3,2i+H  VAPOR, 
6, Fll. 3,Fl0.3,F11.3 
ic»   SOLID-LlCiUlO  T, 
T,1  RHu,  MOL/L 

F8.3,2'.H  SOLID, 
o, Fll. 3, Flu. 3, Fll. 3 
F4.3,2-*H   rilx,    1  IS 
6, Fll. 3, Fi;. 3, Fll. 3 
Fd.Z^i'^^  LIUUID, 
6, Fll. 3, FIG. 3, Fll. 3 
Fi.3,2'*H  SOLID, 
o,Fll.3,FlC.3,Fll.3 
F8.3,2UH    MIX,    Q  IS 

c,  Fll. 3, FlC. 3, Fll. 3 
F'i.3,2'tH  LIQUID, 
o,Fll.3,Fi:.3,Fll.3 
;2H   LIOUIO-VAPOR  T, 
TM  RhC,  MOL/L 
Fi.J,2^H  LIUUID, 
o,F11.3,Fl3.3,F11.3 
Fc.3,2'-H  VAPOR, 
c,Fli.3,rlL.3,F11.3 


,  USL 

K 

H,  J/MQL 
,  /  ) 

FRACTIONAL 
,/) 

♦/^ 
,/  ) 

FRACTIONAL 
♦  /  ) 

,/) 

K 

H,  J/MOL 
,  /) 

FRACTIONAL 
,  /  ) 

,  /  ) 

,  /  ) 

FRACTIONAL 
,/) 

,  /  ) 

K      FROM  TH: 
-\  f  J/MOl 

,  /  ) 

,  /  ) 


a  TP 

5,    J/MOL-K      U,  J/MOL*,//) 
,F3.1,5X 


=  , F3. 1  ,  5X 
, F3. 1 , 5X 
, F3. 1 ,F5  . 1 

=  ,F3.1,F5.1 
,F3. 1 ,F5  .1 


S,  J/MOL-K 
,  F3 . 1 , F5  .  1 

=  ,  F3 . 1  ,F5  .  1 
, F3. 1  ,F5  . 1 
,F3.1,5X 

=  ,F3.i,5x 

,  F  3  .  1  ,  5  X 

,   32   TERM  EG. 
S,  J/MOL-K 
,  F  3  .  1  ,  5  X 

,  F3 . 1 , 5x 


a  TP 

,  J/MOL*,//) 


Q  TP 
,  J/MOL',//) 
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APPENDIX  D.      TEST  ?ROGR?uM  OXYGEN 
3RCGRAM   VALUESIINPUT, OUTPUT) 

THIS  PROGRAM   TESTS   THE   SLUSH   OXYGEN  SUBROUTINE 
CALL  3ATA02 

C0M3INATI0NS   OF   SOLID    -  VAPOR 
PRINT  1 

T=5a . 

Q=C. 

CALL  SL02 (T ,Q, Q. ,PSL , DSL , HSL , SSL , USD 
PRINT  2,T,Q,  PSL»OSL,HSL,SSL,USL 
Q=0.5 

CALL  SL02(T,Q,D.,PSL,[;SL,HSL,SSL,USL) 
PRINT   3,T,Qt  PSL,DSL,HSL,5SL,USL 
Q=l. 

CALL  SL02(T,Q,Q.,PSL,0SL,HSL,SSL,USL} 
PRINT  i*,T,Q,        PSL»OSL,hSL  ,SSL,USL 

To=l. 
Q=C  . 

CALL  SL02(T,Q,TP,PSL,0SL,HSL,SSL,USL) 
PRINT  5,T,Q,TD,PSL,CSLtHSL,SSL,USL 
Q=0.= 

CALL  SL02f T»Q,TP,PSL»OSL,HSL,SSL,USL) 
PRINT  6,T,Q,TP,PSL,DSL,HSL,SSL,USL 
Q=l. 

CALL  SL02(T,Q,TP,PSL,DSL,HSL,SSL,USL) 
PRINT  7,T,Q,TP,PSL,DSL,hSL,SSL,USL 

COMBINATIONS   OF   SOLID   -  LIQUID 

PRINT  8 
T  =  5i».359 
TP=]  , 
q=C. 

CALL  SL02fT,Q,TO,PSL,CSL,HSL,SSL,USL) 
PRINT   9,T,Q,TP ,?SL»DSL,HSL  ,3SL,USL 
Q=u  .5 

CALL  SL02<T,Q,TD,PSL,DSL,HSL,SSL,USL) 

PRINT10,T,Q,TP,PSL,DSL,HSL,SSL,USL 

0=1. 

CALL   SL02(T,0,T3,3SL,DSL  ,HSL,SSL,USL) 
°RINT11,T,Q,TP,PSL,DSL,HSL,SSL,JSL 
T  =  56. 
Q=Q  . 

CALL  SLO2tT,O»0.,PSL,DSL,-^SL,SSL,USL5 
P'INT12, T,3,  3SL,CSL,HSL,SSL,USL 
0=0  .  5 

CALL  SL02{T,0,G.,3SL,OSL,-iSL,SSL»U3L) 
PRINT13,T,0,        3SL  ,  DSL  ,-S;.  ,SSL  ,USL 
Q  =  l. 

CALL    SL0  2(T,::,G.  ,  -  S  L  ,  3  S  L  , --S  L  ,  S  SL  ,  USL ) 
PRINTIU, T ,0,        P3L, DSL,hSl  ,33L,USL 


C0H3INATI0NS   CF  LIQUID 

PRINT  15 
T=5«4,359 


-'J^PZ^   "Rj-^   TME   32   TE  =  M 
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Q  =  0. 

PSL=VFNIT) 

QSL  =  FINOM(PSL»T,i+l.) 
HSL=ENTHflL(PSL,OSL, T) 
USL=HSL-PSL/OSL* 10 1.325 
SSL=rENTROPeOSL,T) 
P1RINT16,T,Q,  PSL,DSL,HSL,SSL 
(3  =  1. 

OSL=FINOM4PSL»T,0. 000  3) 

HSL  =  ENTHAL8PSL  ,OSL »T) 

USL=HSL-PSL/OSL*10i. 32  5 

SSL=ENTROPCOSL,T) 

PRINT17,T,Q,  PSL,OSL,HSLtSSL 

T  =  56. 

Q=G. 

PSL=VPN(T) 

OSL  =  FINDM(PSL«T,i*l.  ) 

HSL=ENTHALCPSL,OSL , T) 

USL=HSL-PSL/OSL*101. 3  25 

S3L=ENTROP60SL,T) 

PRINT16»T,Q,        PSL»OSL ,HSL,SSL 

Q=l. 

OSL=f INOMtPSL.T»0. 0003) 

HSL=eNTHALCPSL»OSL,T) 

USL=HSL-PSL/OSL*101. 3  25 

SSL=ENTROP(IDSL,T» 

PRINT17,T,Q,  PSL,OSL,HSL,SSL 

1  FORMAT (///, 52H1S0LID-VAP0R  T, 
1      ♦  P,    ATM  RHO,  MOL/L 

2  FORMAT  C10X,F8.  3,  2i*H  SOLID, 

1  F13.8,Flif,6,F11.3,F10.3,F11.3 

3  FORMAT  ( 10X»F8.  3  ,2'+H   MIX,    Q  IS 

1  F13,8,Fli».9,F11.3,F10.3,F11.3 
k   FORMAT  (13X,F8.  3, 2i»H  VAPOR, 
1  F13.8  ,Flt+.9,F11.3,F10.3,Fil,3 

5  FORMAT  (10X,F8. 3, 2'+H  SOLID, 

1  F13.8,Fli».6,Fll,3,F10.3,Fll,3 

6  FORMAT  (10X,F8. 3, 2'+H   MIX,    Q  IS 

1  Fl3.8,Fli+.9,Fll,3,F10.3,Fll,3 

7  FORMAT(  1QX,F8.  3,2'+H  VAPOR, 

1  F13,8,Flt,,9,F11.3,FlG.3,F11.3 
3  FORMAT (///, 52H  SOLID-LIQUID  T, 

1      *  P,    ATM  RHO,  MOL/L 

9   FORMAT(10X,F3.3,2'fH  SOLID, 

1  Fl3.8,Flt».6,F11.3,F10.3,Fll,3 
13    FORMAT(10X,F8.3,2i+H   MIX,    Q  IS 

1  F13.8  ,Fltf.6,F11.3,F10.3,Fll,3 

11  F0PMAT(1QX,F8.3,2'+H  LIQUID, 

1  Fl3.8,Fli+.6,Fll,3,F10.3,F11.3 

12  FORMAT  (10X,F3.  3, 2'+H  SOLID, 

1   Fl3.8,Fli».6,Fli,3,rl0.3  ,F11,3 

13  FORMAT  (  10X,Fd  .  3, 2i+H   MIX,   Q  IS 

1  F13.8,Fl<f,6,F11.3,Fl0.3,F11.3 
i<4  •FORMAT(10X,F3.3,2!*H  LIQUID, 
1   F13.8,FH*.6,F11.3,F1C  .3,Fli.3 

15  FORMAT (///,52H   LIQUID-VAPOR  T, 
1     ♦  P,    ATM  RHO,  MOL/L 

16  FORMAT{10X,F3.3,2«»H  LIQUID, 

1  F13.8,F1»+,6,F11.3,F10.3,F11,3 

17  F0PMAT(1JX,F8.3,21+H  VAPOR, 

1  F13.8,Fli+.9,F11.3,Fl0.3,Fll,3 

END 


USL 


,USL 


,USL 


,USL 

K 

H,    J/MOL  S, 
,/) 

FRACTIONAL  =  , 
,/) 

♦ 

,/) 

♦ 

t/) 

FRACTIONAL  =  , 
,/) 

» 

,/) 

K 

H,    J/MOL  S, 

» 

,/) 

FRACTIONAL  =  , 
,/) 

» 

./) 

» 

,/J 

FRACTIONAL  =  , 
,/) 

» 

,/) 

K     FROM   THE  32 
H,    J/MOL  S, 

,/) 
♦  /) 


Q  TP 
J/MOL-K      U,  J/MOL»,//) 
F3.1,5X 


F3.1,5X 
F3.1,5X 
F3.1,F5.i 
F3.1,F5. 1 
F3.1,F5.  1 


J/MOL-K 
F3,1,F5.1 

F3.1,F5.1 

F3.1,F5.1 

F3. 1, 5X 

F3. 1, 5X 

F3.1, 5X 

TERM  EQ. 
J/MOL-K 
F3.1,5X 

F3.1,5X 


Q  TP 
,  J/MOL*,//) 


Q  TP 
,  J/MOL*,//) 
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APPENDIX  E.      CONVERSION  FACTORS,  HYDROGEN 


Temperature 
Pressure 
Specific  Volume 

Internal  Energy,  Enthalpy- 
Entropy,  Specific  Heat 
Thermal  Conductivity 
Viscosity 
Speed  of  Sound 
Molecular  Weight 
Surface  Tension 


1.8  R  =  1  K 

14.695949  psia  =  1  atir.  =  1.01325  x  10^  N/m^ 

0.00794590  ft^/lb     =  1  cm^/mol 
m 

(1  cm^  =  0.001  liter  =  lO"^  m^) 

0.213405  3TU/lb     =  1  J/mol 
m 

0.118558  3TU/lb  R  =  1  J/mol-K 
m 

0.0578176  BTU/ft-hr-R  =  1  mW/cm-K 

0.067196897   lb  /ft-s  =  1  a/cm-s 
m 


3 .2808 


/s  =  1  m/s 
12 


2 . 01594  (on  the  C 
0.5710147  X  10~^ 


=  12.000  scale) 
lb. /in  =  1  dyn/cm 
-5 


:i  dyn  = 


n: 


APPENDIX  F.      CON'VERSION  FACTORS,  OXYGEN 
Temperature  1.3  R  =  1  K 

Pressure  14.695949-  psia  =  1  atm  =  1.  01325  x  10^  N/m^ 

3  3 

Specific  Volume  0.0005005957  ft  /lb     =  1  cm  /g  nol 

m 

Internal  Energy,  Enthalpy  0.013444  6  3TU/lb^  =  1  J/g  nol 

Entropy,   Specific  Heat  0.0074692  3TU/lb  R  =  1  J/q  mol-K 

Thermal  Conductivity  0.0573176  3TU/ft-hr-R  =  1  mW/ cm-K 

Viscosity  0.067196897  lb  /fr-s  =  1  g/cm-s 

Speed  of  Sound  -  3.2808  ft/s  =  1  m/s 

Molecular  Weight  31.99  8  8 

Surface  Tension  0  .  5710147  x  10  ^  lb -/in  =  1  dyn/cm 

(1  dyn  =  10'^  N) 
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